#!/usr/bin/perl
package run_bl2seq;

#  get_alignment_data('/share/tools/blast/bl2seq', $seq1, $seq2
#  take the @results[4,5]
#  expand cDNA using expand_cdna_with_alignment_array
#  e.g. expand_cdna_with_alignment_array($cdna_seq, '-', $results[4]);
#  call once prepare_calculate_yn00
#	call in loop				calculate_yn00
#  call when finished					cleanup_calculate_yn00


require Exporter;


our @ISA = qw(Exporter);
our $VERSION = 1.00;
our @EXPORT = qw(	
					run_bl2seq
					);

use strict;
use warnings;
use IO::File;
use POSIX qw(tmpnam);
use threads;
use threads::shared;
use parse_bl2seq;


my $tmp_file_name1;
my $tmp_file_name2;
my $tmp_file_name3;
my $tmp_file_name4;

sub delete_tmp_files
{
	my $tmpfile_mutex :shared;
	if ($tmp_file_name1)
	{
		lock ($tmpfile_mutex);
		unlink $tmp_file_name1 or die "Couldn't unlink $tmp_file_name1 : $!";
		$tmp_file_name1 = undef;
	}
	if ($tmp_file_name2)
	{
		lock ($tmpfile_mutex);
		unlink $tmp_file_name2 or die "Couldn't unlink $tmp_file_name2 : $!";
		$tmp_file_name2 = undef;
	}
	if ($tmp_file_name3)
	{
		lock ($tmpfile_mutex);
		unlink $tmp_file_name3 or die "Couldn't unlink $tmp_file_name3 : $!";
		$tmp_file_name3 = undef;
	}
	if ($tmp_file_name4)
	{
		lock ($tmpfile_mutex);
		unlink $tmp_file_name4 or die "Couldn't unlink $tmp_file_name4 : $!";
		$tmp_file_name4 = undef;
	}

}

#_________________________________________________________________________________________


#	line_wrap

#			puts a newline every n letters
#			n.b. Do not use for constant parameters eg "$3"

#_________________________________________________________________________________________
sub line_wrap($$)
{
	my @array;
	my $seq = $_[0];
	while (length $seq)
	{
		push @array, substr $seq, 0, $_[1], "";
	}
	return join "\n", @array;
}


# delete files automatically when we exit or die
END
{
	delete_tmp_files();
}	


#_________________________________________________________________________________________

#	get_alignment_data

#		input:

#			$blast_prog			name and path of the blast programme
#			$seq1				sequence in fasta format
#			$seq2				sequence in fasta format

#		returns:

		
#			$score				bitscore
#			$alignment_start	start of the alignment
#			$alignment_length	length of the alignment
#			$e-value			e-value
#			$i_pos				reference to an array of position / length pairs
#			$j_pos				reference to an array of position / length pairs
#			$final_len			length of alignment positions + non_aligned positions
#			$coverage			alignment coverage = $alignment_length / $final_len

#_________________________________________________________________________________________
sub run_bl2seq
{
	my ($blast_prog, $seq1, $seq2, $options) = @_;

	die unless length($seq1);
	die unless length($seq2);
	my $seg = ($options->{'seg'} ? '-F t' : '-F f');
	my $blosom = ($options->{'blosom'} ? "-M BLOSUM$options->{'blosom'}" : '');
	my $e_value = ($options->{'e_value'} ? "-e $options->{'e_value'}" : '');
	my $coils = $options->{'coils'};

	my $len1 = length($seq1);
	my $len2 = length($seq2);
	$seq1 = ">seq1\n".line_wrap($seq1,60)."\n";
	$seq2 = ">seq2\n".line_wrap($seq2,60)."\n";
	
	
	#
	#	prepare sequence files for bl2seq
	#
	{
		my $tmpfile_mutex :shared;
		my ($tmp_fh3, $tmp_fh4);
		lock ($tmpfile_mutex);
		do {$tmp_file_name3 = tmpnam()}
			until $tmp_fh3 = IO::File->new($tmp_file_name3, O_RDWR | O_CREAT | O_EXCL );
		do {$tmp_file_name4 = tmpnam()}
			until $tmp_fh4 = IO::File->new($tmp_file_name4, O_RDWR | O_CREAT | O_EXCL );
		if (!$coils)
		{
			print $tmp_fh3 $seq1 or die "Failed to print\n$@";
			print $tmp_fh4 $seq2 or die "Failed to print\n$@";
		}
	}
	
	# open temp files
	if ($coils)
	{
		my ($tmp_fh1, $tmp_fh2);
		{
			my $tmpfile_mutex :shared;
			lock ($tmpfile_mutex);
			do {$tmp_file_name1 = tmpnam()}
				until $tmp_fh1 = IO::File->new($tmp_file_name1, O_RDWR | O_CREAT | O_EXCL );
			do {$tmp_file_name2 = tmpnam()}
				until $tmp_fh2 = IO::File->new($tmp_file_name2, O_RDWR | O_CREAT | O_EXCL );
		}
		print $tmp_fh1 $seq1 or die "Failed to print\n$@";
		print $tmp_fh2 $seq2 or die "Failed to print\n$@";
		close $tmp_fh1;
		close $tmp_fh2;
		system "/share/tools/coils/ncoils -f  < $tmp_file_name2  1> $tmp_file_name4 2> /dev/null";
		system "/share/tools/coils/ncoils -f  < $tmp_file_name1  1> $tmp_file_name3 2> /dev/null";
	}
	
	# blast away
	my @results = split/\n/, `$blast_prog -i $tmp_file_name3 -j $tmp_file_name4 $blosom -X 50 -p blastp $seg $e_value`;
	
	# can remove temp files now
	delete_tmp_files();

	# analyse blast results
	return (parse_blast_results(@results, $len1, $len2));

}


1;

